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the  MC  results  in  several  aspects  at  low  concentrations.  The  density  functional  theories 
are  successful  in  predicting  the  extremum  of  the  diffuse- layer  potential  as  a  function 
of  surface  charge  present  in  the  simulations. 
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ABSTRACT 

This  paper  reports  density  profiles  and  mean  electrostatic  potentials  of  a  restricted 
primitive  model  double  layer.  Results  for  an  asymmetric  2:1  salt,  predicted  by  three  differ¬ 
ent  nonlocal  free  energy  density  fimctional  approximations,  axe  compared  with  Monte  Carlo 
results  and  with  those  of  Gouy-Chapman  theory.  Results  for  the  diffuse  layer  potential 
axe  also  compared  with  those  of  some  recent  theories  of  the  double  layer.  The  hard-sphere 
contribution  to  the  free-energy  functional  is  based  on  a  nonlocal  generic  model  functional 
proposed  by  Percus.  We  choose  the  Camahan-Starling  equation  of  state  to  calculate  the 
free  energy  of  the  homogeneous  hard-sphere  mixture  which  enters  in  the  hard-sphere  func¬ 
tional.  The  mean  spherical  approximation  for  a  neutral  bulk  electrolyte  is  used  to  model 
the  electrostatic  part  of  the  non-uniform  ion-ion  correlations  present  in  the  interface.  For 
singly  charged  counterions  the  agreement  between  the  density  functional  approximations 
applied  here  and  the  Monte  Carlo  data  is  excellent.  For  doubly  cheirged  counterions  the 
agreement  is  very  good  at  high  concentrations,  but  differs  from  the  MC  results  in  several 
aspects  at  low  concentrations.  The  density  functional  theories  eire  succesfull  in  predicting 
the  extremum  of  the  diffuse-layer  potentied  as  a  function  of  surface  charge  present  in  the 
simulations. 
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I.  INTRODUCTION 


Recently,  a  nonlocal  density  functional  theory  for  the  electrical  double  layer  was 
presented  [1].  A  double  layer  is  the  separation  of  charge  which  occurs  in  an  electrolyte 
in  the  presence  of  an  electrode.  In  these  studies  the  restricted  primitive  model  of  an 
electrolyte  (RPM)  was  used.  In  this  very  simplified  model,  the  ions  are  modeled  as  charged 
hard  spheres  of  equal  diameter  d  and  the  solvent  is  represented  as  a  uniform  medium  with 
dielectric  constant  . 

In  Ref.(l),  the  general  formalism  due  to  Mermin  [2]  is  used  to  construct  a  grand 
potential  functional  for  the  interface.  The  haird-sphere  contribution  to  the  free  energy 
functional  is  based  on  a  nonlocal  generic  model  fimctional  proposed  by  Percus  [3].  It  was 
shown  by  Vanderlick  et  al.  [4],  that  this  generic  functional  can  be  used  as  a  generating 
functional  of  several  known  nonlocal  approaches.  In  particular,  the  generalized  van  der 
Walls  (GVDW)  model  [5],  the  generalized  haxd-rod  model  (GHRM)  [6,7]  and  the  smoothed- 
density  approximation  (SDA)  due  to  Teu’azona  [8]  can  be  generated  with  this  functional. 
This  generic  functional  was  extended  to  mixtures  by  Vanderlick  et  al.  [9]. 

Using  the  GHRM,  the  density  functional  theory  was  applied  in  Ref.(l)  to  a  planar 
interface  for  symmetric  1:1  and  2:2  RPM  electrolytes.  More  recently  [10],  the  same  the¬ 
ory,  combined  with  the  GVDW  model  and  the  SDA  was  also  applied  to  symmetrical  1:1 
and  2:2  salts.  In  all  cases,  use  was  made  of  the  analytic  solution  of  the  mean  spherical 
approximation  MSA  [11]  for  the  bulk  electrolyte  to  approximate  the  electrostatic  part  of 
the  nonuniform  Omstein-Zemike  direct  correlation  function  wich  appears  in  the  theory. 
In  both  works  [1,10],  the  free  energy  of  the  haxd-sphere  system,  which  appears  in  the 
functional,  weis  approximated  with  the  Camahan-Starling  expression  [12].  In  Ref.(lO),  the 
Clausius  equation  of  state  was  also  used. 

The  results  of  Refs. (1,10)  can  be  summarized  as  follows:  the  comparison  of  the 
density  functional  theory  with  the  computer  simulations  of  Torrie  and  Valleau  [13,14]  for 
symmetric  electrolytes  shows  that,  when  combined  with  the  Camahan-Starling  equation 
of  state,  the  theory  is  very  succesful  in  predicting  density  profiles  and  mean  electrostatic 
potential  functions  for  the  RPM  electrolyte.  The  theory  correctly  predicts  the  layering 
effects  present  in  the  computer  simulations  of  1:1  electrolytes  at  1  M  and  high  electrode 
charges  [13],  and  the  charge  reverszJ  phenomena  which  occurs  in  2:2  electrolytes  [14].  The 
three  models  studied  are  in  very  close  agreement  at  low  concentrations  and  low  surface 
charge  densities  but  the  differences  between  them  increases  when  ch2Lrge  or  concentration 
are  increased.  The  comparison  with  the  MC  data  [13]  for  the  diffuse  layer  potentied  eis  a 
function  of  the  surface  charge  density  [10]  shows  that  the  GHRM  is  superior  to  the  other 
two  models  at  high  surface  density  cheirges  for  highly  concentrated  1:1  electrolytes.  This  is 
probably  due  to  a  fortuitous  cancelation  of  errors  in  the  GHRM  since  the  SDA,  which  is  a 
better  model  for  the  non-uniform  hard-sphere  system  [4],  gives  a  better  representation  of 
the  layering  effects  present  at  high  density  charges.  The  three  models  compared  are  only 
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in  qualitative  agreement  with  the  MC  data  for  the  diffuse  layer  potential  as  a  function  of 
surface  charge  density  for  2:2  salts  [10]. 

In  this  ^u:ticle  Vre  report  results  of  the  density  functional  theory  for  2:1  electrolytes  for 
both  positively  and  negatively  charged  electrodes.  As  in  Ref.(lO),  we  present  a  comparison 
between  the  GHRM,  the  GVDW  model  and  the  SDA.  We  also  compare  our  results  with 
those  of  the  Gouy-Chapmaui  [15,16]  theory  as  modified  by  Stem  [17]  (MGC),  and  with  those 
of  some  recent  theories  of  the  double  layer.  These  include  the  most  refined  version  of  the 
modified  Poisson- Boltzmann  (MPB5)  theory  [18],  the  hypemetted  chain/meeui  sphericjil 
approximation  (HNC/MSA)  theory  [19]  and  the  Bom-Green- Yvon  (YBG)  with  the  c'osure 
introduced  by  Caccamo  et  al.  [20]. 

As  pointed  out  previously  by  several  authors  [14,19],  results  for  2:1  salts  are  inter¬ 
esting  not  only  because  their  practical  applications,  but  also  because  the  MGC  theory 
exhibit  both  qualitative  agreement  or  disagreement,  compared  with  the  MC  results,  de¬ 
pending  upon  whether  the  counterion  is  monovalent  or  divzdent. 

This  aU'ticle  is  as  follows.  In  Sec.  II  we  briefly  review  the  theory.  We  report  our 
results  in  Sec.  Ill  and  we  conclude  in  Sec.  IV  with  a  dissension. 


II.  THEORY 


In  Refs.(l)  and  (10)  the  formalism  due  to  Mermin  [2]  was  adopted  to  construct  a 
model  grand  potentieil  functional,  J2({n}),  for  a  mixture  of  charged  particles  in  an  external 
potential  i;a(r).  This  functions^  of  the  density  profiles  can  be  written  in  the  following  form: 

/?f2({n})  {n* })-!-/?  J  d^rnair)vQ{r) 

+  (na(r)/nj»)  -  J  d^r  [n„(r)  -  n^,] 

1  fX 

dX  dX'ca^{r,r';X'l 

Jo 

(2.1) 

where  0  =  1/ifeT,  T  is  the  absolute  temperature  and  k  is  Boltzmann’s  constant.  In  this 
equation  n({n’})  is  the  grand  free  energy  of  a  reference  or  initial  state  at  the  same  tem¬ 
perature  T,  chemical  potentials,  and  volume  V  than  the  interface.  The  reference  state 
was  chosen  to  be  absent  of  any  external  potential.  It  is  a  state  of  uniform  density  in  which 

nj,(r)  =  Ua  =  constant.  (2.2) 

In  the  obtention  of  Eq.(2.1),  a  linear  path  was  chosen  for  the  functional  integration.  The 
path  is  characterized  by  a  single  parauneter  A  which  varies  from  0  to  the  unity.  The  function 


+  X)// [”a(r)-n;»]  [n^(r)-ny  ^ 
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Caffir,  r';  A)  represents  the  Omstein-Zemike  direct  correlation  function  of  a  nonuniform 
mixture  at  a  psurticular  point  on  the  integration  path. 


The  equilibrium  density  profiles,  ne»(r),  minimize  the  greuid  potential  functional, 
n({n}),  for  fixed  Va(r)  and  intermolecular  potentials  Uo^(r,  r').  For  a  restricted  primitive 
model  electrolyte  the  pair  potential  can  be  sepairated  into  its  Coulombic  zind  short-range 
repulsive  contributions 

UQ/j(r,  r')  =  -  r'|)  -f-  u’'(|r  -  r'i),  (2.3) 

where  z^  is  the  valence  of  species  a,  c  is  the  magnitude  of  the  electronic  charge, 

u(r)  =  1/er,  (2.4) 


and 

“"(r)  =  rKd, 

=  0,  r>d.  ^  ^ 

It  is  also  convenient  to  separate  the  Coulombic  aind  short  range  contributions  to  the  ex¬ 
ternal  potentials.  The  impenetrable  wall  produces  a  repulsive  potential  of  the  form, 


u''(i)  =  oo,  X  <  d/2, 

=  0,  X  >  d/2, 


(2.6) 


where  x  is  the  shortest  distance  to  the  wall.  The  surface  charge  density,  cr,  gives  rise  to  a 
one  body  Coulombic  potential  of  the  following  form 

u(x)  =  — 27rcr|x|/e -I- C,  (2.7) 

where  C  is  a  constant  which  determines  the  point  of  zero  potential.  By  functional  difer- 
entiation  of  the  grand  potential  functioned  w.r.t.  n(r),  we  obtain  the  following  expression 
for  the  equilibrium  density  profiles; 

ln{na  (r)/no)  =  -/3  (t;;(r)  -I-  ezarl’ir)) 

+  E  /  -  ""I  /"  (2  8) 

^AF"^({n}) 

6na{r) 


In  the  leist  equation,  rl>{r)  is  the  mean  electrostatic  potential,  Aca/j(r,  r')  is  the  resid- 
ued  peirt  of  the  direct  correlation  function  defined  by 


Co/j(r,  r')  =  -^ZaZ0e^u*^{\r  -  r'|)  +  Co/(r,  r')  -|-  Acasir,  r'). 


(2.9) 
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and  r')  is  the  hard  sphere  contribution  to  the  function  Ca^(r,  r').  The  last  term  on 

the  right  hand  side  of  Eq.(2.8)  is  the  functional  derivative  of  the  excess  free  energy  change 
between  the  reference  state  and  and  the  final  state  due  to  the  hard  sphere  interaction  [1]. 
The  mean  electrostatic  potential  is  a  solution  of  Poisson’s  equation.  In  the  case  of  a  planar 
interface,  Poisson’s  equation  can  be  integrated  to  give 

47rfi  foo 

rp(x)  =  —  dx' (x  -  x')^Zana{x'),  (2.10) 

where  the  condition  of  overall  electroneutrality, 


!  J  dx'  ^  Zariaix')  =  -a, 


(2.11) 


was  imposed. 

In  order  to  obtain  an  expression  for  the  hard-sphere  excess  free  energy  functional 
^^^■^({n}),  we  make  use  of  the  generic  functional  proposed  by  Percus  [3]  as  the  three 
dimensional  generalization  of  the  Helmholtz  free  energy  of  an  inhomogeneous  hard  rod 
system.  This  functional  c«in  be  written  for  mixtures  [9]  as 

F"“"({n})  =  E .  (212) 

where  nQ(r)  and  n^{T)  axe  coarse  grain  densities  which  are  defined  as  spatial  averages  of 
the  local  densities  over  certain  smedl  domains.  For  a  mixture,  the  coarse  grain  densities 
are  defined  as 

”a(r)  =  J  Uc,iT-r'-,{n})n^{r'),  (2.13) 

"a(r)  =  y  dV'ra(r- r';{n})no(r').  (2.14) 

The  weigthting  functions  i/q{t  -  r';  {n})  and  Tq{t  -  r';  {n}),  which  in  the  most  general 
ceise  axe  functionals  of  the  density  distributions  must  be  normalized,  i.e., 

j  d^r'  i/a{r  -  r';  {n})  =  J  d^r'  Ta{T  -  r';  {n})  =  1.  (2.15) 

By  changing  the  specific  form  of  weighting  functions  in  Eqs.(2.13)  and  (2.14)  different 
model  density  functionals  can  be  generated  [4].  In  this  work  we  solved  Ekj.(2.8)  for  the 
three  model  functionals  mentioned  before;  the  GVDW,  the  GHRM,  and  the  SDA  due 
to  Tarazona.  Here  we  briefly  summarize  the  specific  forms  of  the  weigthting  functions 
Ua(r  —  r';  {n})  and  Ta(r  —  r';  {n})  for  these  three  models.  For  the  GVDW  model,  we  have 

i^a(r-r';{n})  =  <5(r-r'),  (2.16) 
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Tair  -  r';  {n})  =  H(d  -  |r  -  r'|)/(47rd^/3),  (2.17) 

where  ^(r)  the  Dirac  delta  function  and  H{r)  is  the  Heaviside  step  function  defined  by 


H(r)  =  1,  r  >  0, 

=  0,  r  <  0. 


(2.18) 


For  the  GHRM  [6,7], 

^,(r  -  r'i  {n})  =  ((d/2  -  |r  -  r’|)/(4i(rf/2)“),  (2.19) 

r„(r  -  r'i  {n})  =  H(d/2  -  |r  -  r'|)/(4T(<i/2)V3),  (2.20) 

whereas  for  the  SDA  [8],  we  have 

i^a(r  -  r';  {n})  =  ^(r  -  r'),  (2.21) 

ra(r  -  r';  {n})  =  u;^°^(|r  -  r'])  +  -  r'|)n5(r) 

^  (2.22) 

0y 

In  the  case  of  a  restricted  primitive  model  electrolyte,  all  the  ions  have  the  same  diameter. 
As  a  result  of  this  simplification,  the  ty’s  become  to  be  independent  on  the  species,  and  we 
can  rewrite  Eq.(2.22)  as 

Ta{r  -  r';  {n})  =  u;(°)(|r  -  r'|)  +  -  r'|)  ^  n^(r) 

^  ^  (2-23) 

0 


Thus,  in  the  SDA,  the  density  dependent  weighting  functions  ro(r— r';  {n})  Eire  expressed  eis 
a  truncated  series  in  the  nonlocal  densities  nj(r).  The  coefficients  ty(°)(|r—r'|),  ly^^^dr— r'|) 
aind  ly^^^dr  —  r'i)  aie  density  independent  weighting  functions  which  satisfy  the  following 
properties: 

J  d^rw^'\\r\)  =  Soi,  »•  =  0,1,2.  (2.24) 

With  Ta(r  —  r';  {n})  defined  in  this  way,  and  using  Eq.(2.14),  n^(r)  can  now  be  expressed 

EIS 

”a(r)  =  ^(r)  +  (--S) 

0  0 


where 


”a^(r)  =  j 


d^r'w^ 


)(|r-r'|)n„(r'),  i  =  0,1,2. 


(2.26) 
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In  the  SDA  the  coefficients  w  are  determined  by  requiring  the  direct  correlation  function 
c^'^(|r  —  r'l)  of  a  uniform  hard-sphere  fluid  obtained  from  the  second  order  functional 
derivative  of  the  excess  free-energy  functional,  ^  be  close  to  that  of  the  Percus- 

Yevick  approximation  for  the  hard-sphere  uniform  fluid  to  the  second  order  in  density  [8]. 
The  prescription  for  the  density  independent  weighting  functions  is 


4Tr(P' 

=  0, 


r  <  d, 
r  >  d, 


w^^\r)  =  0.475  -  0.648  -f  0.113 

=  0.288  -  0.924  -f-  0.764  -  0.187 


(2.27a) 


r  <  d, 

d<r  <  2d,  (2.276) 

r  >  2d, 


«-'2(5)+5Q’  , 


r  <  d, 
r  >  d. 


(2.27c) 


In  the  SDA,  the  total  coarse-grained  density,  n^,  required  in  our  computations  is  given 
by  the  following  root  of  the  quadratic  equation,  Eq.(2.25), 

„  ^  1  -  E,  '‘(So  (£,«!.>))]* 

^  - 2^(7) - ■ 

(2.28) 


According  to  Eq.(2.8),  our  description  of  the  electrical  interface  is  not  complete  with¬ 
out  the  kwoledge  of  the  inhomogeneous  residual  correlation  function  Acaffir,  r';  A)  and  the 
excess  free  energy  per  particle  Fq  of  a  homogeneous  hard-sphere  fluid.  We  approximate 
the  function  AcQ^(r,  r';  A)  with  the  fimction  Aco;j(|r  —  r'|)  of  a  neutral  bulk  electrolyte 
and  use  the  MSA  expression  for  this  function  given  by  [11] 


r  <  d, 
r  >  d. 


(2.29) 


where 


B  =  Kq  -t- 1  —  (1  -f  , 


and  Kp  =  Kod.  The  quantity  kq  is  the  iixverse  Debye  screening  length  given  by 


(2.30) 
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In  this  work  we 
per  particle  of  a 
the  expression 


use  the  quasiexact  Camahan-Starling  formula  for 
homogeneoxis  heird-sphere  fluid  [12].  For  the  RPM 


m{n) 


t){4  -  3t)) 
(1  —  eta)'^  ’ 


T]  =  nnd^/6, 


the  excess  free  energy 
electrolyte  we  can  use 

(2.31) 


for  the  one  component  haxd-sphere  fluid. 


III.  RESULTS 


In  this  work  we  have  carried  out  calculations  for  2:1  electrolytes  at  bulk  concentrations 
of  0.005,  0.05  and  0.5  M.  For  the  two  latter  concentrations  we  considered  both  positive  and 
negative  surface  charge  densities  <7.  For  concentration  of  0.005  M  the  MC  data  reported 
[14]  are  for  negative  values  of  a  only.  Since  we  are  interested  in  a  compjirison  of  our 
results  with  the  data  obtained  via  computer  simulation,  we  restricted  our  cailculations  at 
bulk  concentration  of  0.005  M  to  negative  values  of  a.  For  positive  surface  charge,  the 
counterions  zu'e  singly  charged  ions  ,  while  for  negative  surfeice  charge  the  counterions  are 
doubly  charged. 

In  order  to  solve  Exi.(2.6)  we  employed  the  finite-element  numerical  techniques  de¬ 
scribed  elsewhere  [21].  As  in  Refs.(l)  and  (10),  we  used  quadratic  Lagrange  interpolating 
polynomials  as  basis  functions  for  solving  Eq.(2.6)  for  the  GHRM  and  linear  or  "chapeau" 
functions  for  solving  the  same  equation  for  the  GVDW  model  and  the  SDA.  We  generated 
solutions  in  the  domain  d/2  <  x  <  R,  where  z  is  the  disteince  to  the  wall  and  R  is  the 
cutoff  distance  used  as  the  upper  limit  in  the  integrals.  A  uniform  mesh  of  N  points  was 
used  in  2dl  the  cases.  Both  N  and  R  are  strongly  dependent  on  the  concentration  of  the 
bulk  electrolyte. 


A  Newton  iterative  scheme  was  used  to  solve  the  set  of  non-linear  algebraic  equations 
for  the  values  of  the  reduced  density  profiles  at  the  positions  of  the  nodes:  gai  =  no(z,  )/no. 
The  iterative  process  is  continued  until  the  Euclidean  norm  of  the  updates  after  the  itera¬ 
tion  k-f-1  becomes  less  than  10“^°,  i.e.. 


E.  {si'*"  -  si") 

2N 


2  1 


i 

<  io~^°. 


(3.1) 


As  a  second  internal  test  of  our  calculations  we  also  have  cheked  the  agreement  of  the 
electroneutrality  condition,  Eq.(2.11).  This  equation  was  satisfied  to  at  least  five  significant 
figures  in  all  our  calculations. 


We  found  convenient  to  use  dimensionless  parameters.  The  dimensionless  surface 
chMge  density  is  a*  =  adP /t.  The  dimensionless  mean  electrostatic  potential  is  i/’*(x)  = 
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0eip{x).  The  distance  x  is  reduced  by  the  diameter  d.  As  in  our  previous  work,  in  order  to 
compare  with  the  MC  data  of  Valleau  and  co-workers  [14],  we  fixed  the  plasma  parameter 
to  r*  =  ^e^fed  =  1.6809.  This  value  of  F*  corresponds  to  r=298  K,  c  =  78.5  and  d= 
4.25A. 

Our  results  for  the  diffuse  layer  potential  at  several  values  of  the  surface  charge,  cr*  , 
for  2:1  electrolytes  are  summarized  in  Table  I.  For  comparison,  in  the  same  table  we  also 
display  the  MC  results  of  Torrie  and  VaJleau  [14],  the  results  of  the  MGC  theory  [15-17], 
the  BGY  theory  of  Caccamo  et  al.  [20],  those  corresponding  to  the  MPB5  theory  [18], 
and  HNC/MSA  results  due  to  Lozada-Cassou  zmd  Henderson  [19].  For  positive  values 
of  <7*,  the  agreement  between  the  results  of  the  density  functional  theories  presented  in 
this  article  and  the  MC  data  is  excellent.  The  GHRM  theory  is  slightly  superior  to  the 
other  two  density  functional  theories  for  positive  a*  .  For  negative  siurface  charges,  the 
agreement  between  our  results  and  the  MC  data  is  only  fair  at  low  concentrations  but 
considerable  improves  at  the  higher  concentration  of  0.5  M.  The  GVDW  theory  and  the 
SDA  Eire  slightly  better  than  the  GHRM  for  the  negative  values  of  a*  investigated  in  this 
work. 


A  comparison  of  the  results  of  the  three  density  functional  theories  with  the  MGC 
theory  and  MC  data  is  presented  in  Fig.  1.  This  figure  shows  the  dimensionless  diffuse 
layer  potential,  ^”*(0)  =  0eiJ^(O)  as  a  function  of  the  electrode  charge  density  a*  for  c  = 
0.005,  0.05  and  0.5  M.  The  MGC  theory,  which  neglects  the  effects  due  to  the  finite  size  of 
the  ions,  overestimates  the  thickness  of  the  double  layer  and  predicts  a  monotonic  rise  for 
the  diffuse  layer  potential  as  a  function  of  <r*.  As  has  been  pointed  out  by  other  authors 
[14],  for  positive  surface  chEirge,  the  behavior  of  the  diffuse  layer  potentiaJ  as  a  function 
of  <7*  is  very  reminiscent  of  that  of  1:1  electrolytes.  For  negative  surface  charges,  the 
monotonic  behavior  of  the  MGC  theory  is  in  significant  disagreement  with  the  MC  data, 
since  the  simulation  results  present  a  minimum  in  the  diffuse  layer  potential  as  a  function 
of  cr*.  The  density  functional  theories  predict  the  presence  of  that  minimun.  Since  the 
results  of  the  GVDW  theory  are  undistinguishable  from  those  of  the  SDA  with  the  scales 
used  in  all  our  plots,  in  what  follows  we  only  show  one  line  for  these  two  approximations. 
For  c  =  0.005  M,  the  results  of  the  three  density  functional  theories  are  undistinguishable 
on  the  plot. 

We  proceed  now  to  examine  the  information  contained  in  the  ionic  density  profiles  and 
mean  electrostatic  potential  profiles.  For  positive  surface  densities  and  low  concentrations, 
the  behavior  of  the  double  layer  is  very  similau-  «,o  that  found  in  1:1  systems.  Thus  we 
describe  the  behavior  at  higher  concentrations  where  more  interesting  phenomenae  are 
present.  In  Fig.  2  we  display  the  density  profiles  at  0.5  M  and  a*  =0.20.  The  solid  lines 
correspond  to  the  GHRM  theory,  while  the  dot-dashed  lines  correspond  to  the  GVDW 
theory  zmd  SDA.  The  circles  correspond  to  the  MC  data  of  Torrie  and  Valleau  [14].  The 
dashed  lines  Eire  the  MGC  density  profiles  calculated  with  the  Euialytic  solution  for  2:1 
RPM  electrolytes  found  by  Grahame  [22].  As  shown  in  Fig.  2,  the  three  density  functional 
theories  predict  a  counterion  profile  with  a  very  shallow  minimum  axoimd  x  =  2.5d.  To 
this  minimum,  corresponds  a  mEiximum  in  the  coion  density  profile.  Note  that  the  position 
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of  the  plate  has  been  shifted  to  a:  =  —d/2  in  all  our  plots.  The  MC  results  seem  to  present 
the  same  oscillatory  behavior,  but  it  is  partially  hiden  by  the  statistical  noise.  As  usual, 
the  classical  MGC  theory  predicts  monotonic  density  profiles.  The  oscillatory  behavior  is 
zilso  apparent  in  the 'mean  electrostatic  potential  as  a  function  of  the  distance  to  the  waJl. 
In  Fig.  3  we  can  see  how  the  MC  data  for  the  fimction  0(i)  have  a  change  in  sign  around 
X  =  1.5d  and  present  a  shallow  minimun.  The  three  density  functional  theories  studied 
predict  quite  well  the  location  and  size  of  this  minimum  in 

We  now  shift  our  attention  towards  systems  with  negative  surface  charges.  This 
situation  corresponds  to  a  system  with  doubly  charged  counterions  and  singly  charged 
coions.  For  systems  with  low  surface  density,  the  behavior  of  the  density  profiles  and 
potential  profile  is  monotonic.  This  is  the  case  for  c  =  0.005  M  and  o*—  -0.05,  see  Fig. 
4,  where  the  results  of  the  three  density  functional  theories  for  the  fimction  ^(x)  are 
coincident  up  to  the  scale  presented.  Even  when  the  results  obtained  with  the  density 
functional  theories  are  in  considerable  improvement  over  the  MGC  theory,  appreciably 
overestimate  the  interfacial  thickness. 

At  higher  concentration  and  higher  surface  ch2irges  the  MC  profiles  are  no  longer 
monotoD’C.  In  Fig.  5  we  present  MC  density  profiles  at  c  =  0.05  M  and  <7*  =  —0.20. 
The  coion  density  profile  obtained  in  the  MC  simulation  [14]  has  a  peak  around  x  =  1.5d. 
In  the  same  figure  we  display  the  profiles  obtained  with  the  density  functional  theories 
studied  in  this  paper  and  those  of  the  MGC  theory  The  profiles  predicted  by  the  density 
functional  theories  are  still  monotonic  but  show  a  clear  tendency  towards  an  oscillator^' 
behavior.  As  before,  the  differences  between  the  GVDW  emd  the  SDA  are  negligible. 
Differences  between  the  GHRM  theory  and  the  other  two  density  functional  theories  are 
also  very  small  in  this  case.  The  three  theories  produce  coion  density  profiles  which  are 
undistinguishable  in  the  figure.  At  the  larger  surface  charge  of  a*  =  —0.284  eind  the  same 
concentration  c  =  0.05  M,  the  MC  data  show  a  very  small  charge  inversion  in  between 
X  =  2.5<i  and  4.5<i;  see  Fig.  6.  The  profiles  predicted  by  the  density  functional  theories  are 
still  monotonic  and  the  differences  between  them  continue  to  be  very  small.  Of  course, 
the  MGC  theory  profiles  are  monotonic.  The  chzirge  inversion  phenomenae  present  in  the 
MC  data  at  these  conditions,  produces  an  oscillatory  behavior  in  the  mean  electrostatic 
potential  profile;  see  Fig.  7. 

The  oscillatory  behavior  in  the  MC  data  is  more  evident  at  higher  concentrations. 
In  Fig.  8  we  present  a  comparison  of  the  MC  density  profiles  with  those  obtained  with 
the  density  functionad  theories  studied  here  at  c  =  0.5  M  and  cr*  =  —0.1704.  The  charge 
reversal  phenomenae  is  very  pronounced  at  these  conditions.  For  comparison,  in  the  same 
figure  we  display  the  results  of  the  MGC  theory.  The  agreement  between  the  density 
functional  theories  and  the  simulation  results  is  good  here.  Again,  the  differences  between 
the  three  density  functioned  theories  are  very  small.  The  cheurge  reversal  found  at  c  =  0.5  M 
and  <7*  =  —0.1704  produces  a  relatively  deep  minimum  in  the  mean  electrostatic  potential 
profile.  In  Fig.  9  we  compare  the  results  of  the  density  functioned  theories  for  the  function 
V’(x)  with  the  MC  results.  The  agreement  between  the  theories  theories  studied  here  and 
the  simulation  results  is  good,  but  the  theories  clearly  underestimate  the  deepness  of  the 
minimum. 
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IV.  SUMMARY 


In  this  article  we  applied  three  different  nonlocal  free  energy  density  functional  ap¬ 
proximations  to  the  problem  of  a  double  layer  consisting  of  a  2:1  RPM  electrolyte  solution 
in  the  presence  of  a  planar  electrode.  Previously,  we  applied  these  theories  to  a  double 
layer  in  equilibrium  with  symmetrical  1:1  and  2:2  salts  [1,10].  These  nonloced  density 
functional  approximations  are  based  on  the  generic  functional  proposed  by  Percus  [3]  and 
generalized  to  mixtures  by  Vandelick  et  al.  [9].  This  functional  attempts  to  account  for 
the  nonlocal  effects  due  to  the  finite  size  of  fiuid  particles  by  incorporating  coarse-grain 
densities.  The  approximations  studied  here  provide  three  different  models  of  weighting 
functions  for  the  specification  of  the  coarse-grain  densities.  Following  Refs.(l)  and  (10), 
we  use  the  Camahan-Starling  [12]  equation  of  state  to  approximate  the  free  energy  of  a 
homogeneous  hard-sphere  mixture.  This  latter  quantity  is  required  as  an  imput  by  the 
functional.  The  residual  electrostatic  part  of  the  inhomogeneous  direct  correlation  func¬ 
tions  are  approximated  with  those  of  the  neutral  bulk  electrolyte  which  is  in  equilibrium 
with  the  interface.  We  use  the  analytical  solution  of  the  MSA  for  bulk  electrolytes  [ll]  for 
the  latter. 

The  density  functional  approximations  compared  in  this  article,  are  quantitatively 
correct  in  their  predictions  of  density  profiles  and  mean  electrostatic  potential  functions 
of  2:1  electrolytes  when  the  counterions  are  singly  chairged.  The  values  of  the  diffuse  layer 
potential  obtained  with  these  theories  for  positive  cr,  are  in  excellent  agreement  with  the 
MC  data  within  the  intervals  of  composition  and  surface  charge  explored.  The  GHRM  is 
slightly  better  than  the  other  theories  studied  in  this  work  for  singly  charged  counterions. 
This  result  is  reminiscent  of  that  found  for  1:1  electrolytes  [10].  At  very  high  charges, 
one  would  expect  the  appearence  of  the  layering  phenomenon  observed  in  1:1  salts.  The 
differences  between  the  theories  compared  and  the  simulatio  results,  grow  very  slowly 
when  (7  is  increased. 

The  theories  compared  in  this  article  correctly  predict  the  qualitative  asynunetry  with 
respect  to  the  sign  of  the  surfeu:e  charge  found  in  the  MC  data  of  the  diffuse  layer  potential. 
However,  only  at  high  concentration  the  theories  give  a  satisfactory  quantitative  prediction 
of  this  property  when  the  counterions  are  doubly  charged.  At  low  concentrations  and 
small  negative  surface  charges,  the  theories  correctly  predict  the  monotonic  behavior  of  the 
density  aind  mean  electrostatic  potenti2d  profiles,  but  fail  to  predict  the  oscillatory  behavior 
present  in  the  MC  results  when  the  surface  charge  is  increaised.  At  high  concentration  the 
density  and  mean  electrostatic  potential  profiles  obtained  in  this  work  show  the  charge 
reversed  effect  in  agreement  with  the  MC  simulations. 

For  singly  charged  counterions,  the  theories  studied  here  compete  in  accuracy  with  the 
results  obtained  by  Plischke  emd  Henderson  [23]  by  solving  the  inhomogeneous  Omstein- 
Zemike  (OZ)  equation  for  the  pair  correlation  function  together  with  the  Lovett-Mou-Buff- 
Wertheim  equation  and  the  hypernetted  chain  closure.  For  doubly  charged  counterions. 


11 


the  OZ/LMBW  theory  shows  to  be  superior  by  predicting  the  oscillatory  behavior  present 
in  the  density  profiles  and  mean  electrostatic  potential  at  low  concentrations  and  high 
surface  chzirges. 

i 

The  GVDW  model  is  the  zero  order  approximation  of  the  SDA.  Thus,  the  very  small 
differences  found  between  the  predictions  of  these  two  theories  within  the  intervals  of 
composition  and  a  studied,  for  both  positive  and  negative  values  of  cr,  show  that  the  zero 
order  term  in  the  weighting  functions  of  the  SDA  is  extremely  dominant  over  the  other 
two  terms  in  this  regime. 
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‘TABLE  I.  Diffuse  layer  potential  for  2:1  electrolytes 


MC 

MGC 

BGY 

MPB5 

HNC/MSA 

GHRM 

GvdW 

SDA 

0.005 

H 

-0.01 

-0.91(0.02) 

-0.959 

.  . 

.  . 

•  - 

-0.947 

-0.947 

-0.947 

-0.02 

-1.37(0.01) 

-,1.518 

-1.485 

-1.484 

-1.484 

-0.05 

-1.87(0.03) 

-2.381 

-  “ 

-2.233 

-2.228 

-2.228 

0.05 

M 

-0.05 

-1.05(0.01) 

-1.315 

-1.11 

-1.05 

-1.158 

-1.160 

-1.156 

-1.156 

-0.0975 

-1.21(0.06) 
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FIGURE  CAPTIONS 


Fig.  1  Reduced  diffuse  layer  potential,  /3eil>(0),  as  a  function  of  the  charge  density,  a<P /e, 
for  2:1  electrolytes.  Solid  lines  represent  the  results  of  the  GHRM,  dot  dashed  lines 
the  results  of  GVDW  and  SDA.  The  results  displayed  are  for  0.005,  0.05  and  0.5 
M.  Solid  squares,  open  circles,  and  solid  circles  are  the  corresponding  MC  results. 
Dashed  lines  are  the  results  of  the  MGC  theory. 

Fig.  2  Reduced  density  profiles  n(x)/n  for  a  2:1  electrolyte  at  c  =  0.5  M  and  <t*  =  0.2.  The 
dots  are  the  MC  results.  The  dashed  lines  correspond  to  the  MGC  theory,  the  solid 
lines  to  the  GHRM,  and  the  dot  dashed  lines  to  the  GVDW  and  SDA.  Note  that  the 
wall  is  at  X  =  —d/2. 

Fig.  3  Reduced  mean  electrostatic  potential  profile  at  c  =  0.5  M  and  a*  =  0.2  Symbols  as 
in  Fig.2. 

Fig.  4  Reduced  mean  electrostatic  potential  profile  for  a  2:1  RPM  electrolyte  at  c  =  0.005  M 
and  <7*  =  —0.05.  The  dots  are  the  MC  results  of  Ref. (13),  The  dashed  line  correspond 
to  the  MGC  theory,  and  the  solid  line  to  the  density  functional  approximations  used 
in  this  work.  At  the  low  concentration  and  surface  charge  corresponding  to  this  figure, 
the  results  of  the  three  density  functional  approximations  are  undistinguishable  on 
the  plot. 

Fig.  5  Reduced  density  profiles,  n(x)/n,  for  a  2:1  electrolyte  at  c  =  0.05  M  and  a*  =  —0.2. 
The  dots  are  the  MC  results.  The  solid  lines  correspond  to  the  GHRM,  and  the  dot 
dashed  lines  to  the  GVDW  and  SDA.  Dashed  lines  are  the  results  of  the  MGC  theory. 

Fig.  6  Reduced  density  profiles,  n{x)/n,  for  a  2:1  electrolyte  at  c  =  0.05  M  and  a*  =  —0.284. 
All  symbols  as  in  Fig.  5. 

Fig.  7  Reduced  meaui  electrostatic  potential  profile  for  a  2:1  electrolyte  at  c  =  0.05  M  rind 
a*  =  —0.284.  All  symbols  as  in  Fig.  5. 

Fig.  8  Reduced  density  profiles,  n(x)/n,  for  a  2:1  electrolyte  at  c  =  0.5  M  and  a*  =  —0.1704. 
The  dots  are  the  MC  results.  The  solid  lines  correspond  to  the  GHRM,  amd  the  dot 
dashed  lines  to  the  GVDW  and  SDA.  Dashed  lines  are  the  results  of  the  MGC  theory. 
The  density  functionaJ  theories  predict  the  charge  reversal  phenomenon. 

Fig.  9  Reduced  mean  electrostatic  potential  profile  for  a  2:1  electrolyte  at  c  =  0.5  M  and 
a*  =  —0.1704.  The  dots  are  the  MC  results.  The  solid  line  correspond  to  the  GHRM, 
and  the  dot  dashed  line  to  the  GVDW  and  SDA.  The  dashed  line  corresponds  to  the 
MGC  theory. 
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